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ABSTRACT 

We present the discovery by the WASP-South survey, in close collaboration with 
the Euler and TRAPPIST telescopes, of WASP-121 b, a new remarkable short- 
period transiting hot Jupiter, whose planetary nature has been statistically val¬ 
idated by the FASTIS software. The planet has a mass of -^Jupj a 

radius of 1.865 ± 0.044 i?jup, and transits every 1.2749255to ooooo25 ac¬ 
tive F6-type main-sequence star (14=10.4, 1.353lQ;g7g 1-458 ± 0.030 Rq, 

Teff = 6460 ± 140 K). A notable property of WASP-121 b is that its orbital semi¬ 
major axis is only ~1.15 times larger than its Roche limit, which suggests that the 
planet might be close to tidal disruption. Furthermore, its large size and extreme ir¬ 
radiation (^7.110® erg s“^cm“®) make it an excellent target for atmospheric studies 
via secondary eclipse observations. Using the TRAPPIST telescope, we indeed detect 
its emission in the z'-band at better than ^4cr, the measured occultation depth being 
603 ± 130 ppm. Finally, from a measurement of the Rossiter-McLaughlin effect with 
the CORALIE spectrograph, we infer a sky-projected spin-orbit angle of 257.815)5 
This result indicates a significant misalignment between the spin axis of the host star 
and the orbital plane of the planet, the planet being in a nearly polar orbit. Such 
a high misalignment suggests a migration of the planet involving strong dynamical 
events with a third body. 

Key words: planetary systems - stars: individual: WASP-121 - techniques: photo¬ 
metric - techniques: radial velocities - techniques: spectroscopic 
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1 INTRODUCTION 

Most of the transiting exoplanets found by ground-based 
transit surveys (e.g. WASP, Pollacco et al. 2006; HATNet, 
Bakos et al. 2004) are Jovian-type planets with orbital pe¬ 
riods of just a few days, these planets being the easiest to 
detect for such surveys. The orbital period distribution of 
these so-called “hot Jupiters” is not smooth and presents 
a pile-up around periods of ~3-4 days (see e.g. Gumming 
2010). While the long-period drop-off can be explained by 
a lower transit probability for these systems combined to 
a selection effect, the reduced number of planets in orbital 
periods less than 2 days is definitely real, being seen in both 
ground- (e.g. WASP, Hellier et al. 2012) and space-based 
(e.g. Kepler, Howard et al. 2012) transit surveys, as well as 
in radial velocity surveys (see e.g. Marcy et al. 2004). 

Ford & Rasio (2006) suggested that the lower edge of 
the pile-up is defined not by an orbital period, but rather 
by a tidal limit, and found that the inner cutoff is actually 
close to twice the Roche limit (ur)^. This can be naturally 
explained if planets were initially scattered into highly ec¬ 
centric orbits with short pericenter distances from much fur¬ 
ther out, due e.g. to planet-planet interactions (e.g. Rasio 
& Ford 1996, Weidenschilling & Marzari 1996, Moorhead 
& Adams 2005, Chatterjee et al. 2008) and/or Kozai cycles 
(e.g. Kozai 1962, Lidov 1962, Wu & Murray 2003, Fabrycky 
& Tremaine 2007), and later circularized via tidal dissipa¬ 
tion. On the contrary, they argued that this result is incon¬ 
sistent with a disk-driven migration scenario (e.g. Goldreich 
& Tremaine 1980, Lin & Papaloizou 1986, Tanaka et al. 
2002, Lubow & Ida 2010), as the inner edge of the orbital 
period distribution should then be right at the Roche limit. 
The observed distribution of orbital obliquities, with many 
planets found on misaligned or retrograde orbits (e.g. Triaud 
et al. 2010), also supports dynamical migration processes in¬ 
volving a third perturbing body, rather than disk migration. 

The finding of several hot Jupiters with orbital separa¬ 
tions a lower than 2 ur, such as WASP-12 b (a/aR~1.09, 
Hebb et al. 2009), WASP-19 b (a/aR~1.08, Hebb et al. 
2010), WASP-103b (a/oR~1.16, Gillon et al. 2014), 
OGLE-TR-56b (a/aR~1.23, Konacki et al. 2003, Adams 
et al. 2011), and WTS-2b (a/aR~1.27, Birkby et al. 2014), 
challenged the scattering scenario as these planets would 
have been destroyed or completely ejected from their sys¬ 
tems if they had been directly scattered to such short peri¬ 
center distances (see e.g. Guillochon et al. 2011). However, 
Matsumura et al. (2010) showed that these extreme orbits 
can still result from the scattering model, assuming first a 
scattering into an eccentric orbit beyond 2 or, followed by a 
slow inward migration and circularization through tidal dis¬ 
sipation inside the planet mainly until reaching ~2 or, and 
from then tidal decay through tidal dissipation inside the 
star only. The speed of the final tidal decay depends on the 
tidal dissipation efficiency of the star, which is parameter¬ 
ized by Q/, the stellar tidal dissipation factor. Despite being 
an essential parameter in the theory of stellar tides, Q/ is 
still poorly constrained, with estimates based on theoretical 
and observational studies ranging from to 10® to 10® (see 
e.g. Jackson et al. 2008, Ogilvie 2009, Penev et al. 2012). 

^ l.e. the critical orbital separation inside which a planet would 
lose mass via Roche lobe overflow. 


Planets in the a/aR<2 regime are thus key objects to 
further advance our understanding of how tidal forces in¬ 
fluence the orbital evolution of close-in giant planets. Fur¬ 
thermore, these planets being highly irradiated due to their 
proximity to their host stars, they are also generally favor¬ 
able targets for atmospheric studies via secondary eclipse 
observations (see e.g. Seager & Deming 2010, Anderson et al. 
2013, Gillon et al. 2010). They thus provide us with a unique 
opportunity to study the relationship between the observed 
atmospheric properties of hot Jupiters and their tidal evolu¬ 
tion stage. In this paper, we report the discovery of a new hot 
Jupiter of this rare kind by the WASP survey, WASP-121 b, 
which orbits a 10.4 U-magnitude F-type star at just ~1.15 
times its Roche limit. 

Section 2 presents the WASP discovery photometry, as 
well as the follow-up photometric and spectroscopic observa¬ 
tions that we used to confirm and characterize the system. 
In Section 3, we describe the spectroscopic determination 
of the stellar properties and the derivation of the system 
parameters through a combined analysis of our photomet¬ 
ric and spectroscopic data. The statistical validation of the 
planet is then described in Section 4. Finally, we discuss our 
results in Section 5. 

2 OBSERVATIONS 

2.1 WASP transit detection photometry 

The WASP transit survey is operated from two sites with 
one for each hemisphere: the Observatorio del Roque de 
los Muchachos in the Ganary Islands in the North and 
the Sutherland Station of the South African Astronomical 
Observatory (SAAO) in the South. Each facility consists 
of eight Ganon 200mm f/1.8 focal lenses coupled to e2v 
2048x2048 pixels GGDs, which yields a field of view of 
450 deg® for each site with a corresponding pixel scale of 
13.7”/pixel. Further details of the instruments, survey, and 
data reduction procedures can be found in Pollacco et al. 
(2006), while details of the candidate selection process can 
be found in Gollier Gameron et al. (2007). 

The host star WASP-121 (1SWASPJ071024.05- 
390550.5 = 2MASS07102406-3905506, V=10.4, K=9A) 
was observed by the WASP-South station (Hellier et al. 
2011) from 2011 Oct 28 to 2012 Mar 29, leading to the 
collection of 9642 photometric measurements. These data 
were processed and searched for transit signals, as described 
in Gollier Gameron et al. (2006), leading to the detection of 
periodic dimmings of about 1.6% with a period of 1.27 days. 
Fig. 1 presents the WASP photometry folded on the best-fit 
transit ephemeris. 

The sine-wave fitting method described in Maxted 
et al. (2011) was used to search for periodic modulation in 
the WASP photometry of WASP-121 that would be caused 
by the combination of stellar activity and rotation, but 
no periodic signal was found above the mmag amplitude. 
This analysis was performed over the frequency interval 
0-1.5 cycles/day at 8192 evenly spaced frequencies. 

2.2 Spectroscopy and radial velocities 

The GORALIE spectrograph, mounted on the 1.2m Euler- 
Swiss telescope at the ESO La Silla Observatory (Ghile), 
was used to gather eighty-nine spectroscopic measurements 
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Figure 1. WASP photometry for WASP-121 folded on the best- 
fit transit ephemeris from the transit search algorithm presented 
in Collier Cameron et al. (2006), and binned per O.OOSd intervals. 
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Table 1. CORALIE radial-velocity (RV) measurements for 
WASP-121. The uncertainties (crv) S'l'e the formal errors (i.e. 
with no added jitter). The uncertainties on the CCF bisector span 
(BS) and FWHM values are 2.5 o-rv This table is available in its 
entirety via the CDS. 


of WASP-121 between 2013 Sep 11 and 2015 Jan 12 (we 
note that the optical fibre feeding the instrument was re¬ 
placed in Nov 2014). Among these spectra, nineteen were 
obtained dnring the transit of 2014 Dec 24 and eighteen dur¬ 
ing the transit of 2015 Jan 12, with the aim of measuring the 
Rossiter-McLaughlin (RM) effect (Rossiter 1924, McLaugh¬ 
lin 1924). WASP-121 was indeed considered an interesting 
target for such measurements, as its high projected rotation 
velocity u* sin i* of 13.5 ± 0.7 km s“^ (see Section 3.1), 
combined with the observed transit depth, was expected to 
yield a RM effect with a semi-amplitude ~135 ms“^. 

Radial velocities (RVs) were computed from the spec¬ 
tra by weighted cross-correlation (Pepe et al. 2002), using a 
numerical G2-spectral template that provides optimal pre¬ 
cisions for late-F to early-K dwarfs (Table 1). A preliminary 
orbital analysis of the RV time-series revealed a 1.27 days 
periodic variation (see the top panel of Fig. 8), in phase 
with the WASP photometry, and with a semi-amplitude 
~180 ms“^ compatible with a planetary-mass companion 
(see Fig. 2a or b). The RM effect was found to have a sur¬ 
prisingly low total amplitude ~70 ms“'^, suggesting that 
the planetary orbit was likely polar (see Fig. 3). It also ap¬ 
peared that the star exhibits an especially high scatter in its 
RV residuals: the standard deviation of the best-fit residu¬ 
als is 67.1 ms“^ for a circular model and 66.0 ms“^ for an 
eccentric model, while the average RV error is 30.7 ms“^. 
We consider further the origin of this high jitter in Section 
3.2. 

The cross-correlation function (CCF) bisector span 
(Queloz et al. 2001) and FWHM values are plotted in Fig. 
2c and 2d, respectively. Both present large variations, their 
standard deviations being 190.2 ms“^ and 245.9 ms“^, re¬ 
spectively, while their average error (calculated as 2.5 times 
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Figure 2. a) CORALIE radial velocities (RVs) for WASP-121 
phase-folded on the best-fit transit ephemeris, along with the 
best-fit circular model and residuals (jitter is not included in 
the error bars), b) Same as the top panel but here first-order 
polynomial functions of the CCF bisector span and FWHM were 
subtracted from the RVs (see Section 3.3). The scatter in the RV 
residuals is significantly reduced, c) Change in the CCF bisector 
span as a function of orbital phase, d) Change in the CCF FWHM 
as a function of orbital phase. 
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Figure 3. Top: Zoom on the Rossiter-McLaughlin effect observed 
with CORALIE. The RVs obtained during the transit of 2014 Dec 
24 are plotted in blue, while the RVs obtained during the transit 
of 2015 Jan 12 are plotted in red. The superimposed, solid black 
line is our best-fit model. Bottom: Corresponding residuals. 
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Table 2. Summary of the follow-up eclipse photometry obtained for WASP-121. For each lightcurve, this table shows the date of 
acquisition (UT), the used instrument, the eclipse nature, the filter and exposure time, the number of data points, the selected baseline 
function, the standard deviation of the best-fit residuals (unbinned and binned per intervals of 2 min), and the deduced values for /3m, 
/3r and CF = /3m x I3r (see Section 3.3 for details). For the baseline function, p(e^) denotes, respectively, a A-order polynomial function 
of time (e = t), airmass (e = a), PSF full-width at half maximum (e = /), background (e = 6), and x and y positions (e = xy). For the 
TRAPPIST data, the symbol o denotes an offset fixed at the time of the meridian flip. 


Date (UT) 

Instrument 

Eclipse nature 

Filter 

Texp 

(S) 

Np 

Baseline function 

(J 

(%) 

cri20s 

(%) 

dm 

dr 

CF 

2013 Dec 09 

TRAPPIST 

Transit 

Sloan- 2 ' 

13 

763 

p(ti) + 0 

0.26 

0.12 

1.12 

1.52 

1.70 

2013 Dec 25 

TRAPPIST 

Occultation 

Sloan- 2 :' 

13 

902 

p(a^+xy^) + 0 

0.37 

0.17 

1.44 

1.07 

1.54 

2013 Dec 30 

TRAPPIST 

Occultation 

Sloan- 2 :' 

13 

653 

p{P+xy^) + 0 

0.29 

0.14 

1.44 

1.44 

2.06 

2014 Jan 01 

TRAPPIST 

Transit 

Sloan- 2 :' 

13 

765 

p{P+xy^) + 0 

0.25 

0.13 

1.13 

2.35 

2.66 

2014 Jan 13 

TRAPPIST 

Occultation 

Sloan- 2 :' 

13 

867 

p(a^^xy^) -h 0 

0.27 

0.12 

1.07 

1.56 

1.66 

2014 Jan 20 

EulerCam 

Transit 

Gunn-r' 

50 

235 

p{t^+n 

0.10 

0.07 

1.38 

1.28 

1.76 

2014 Jan 24 

EulerCam 

Transit 

Gunn-r' 

50 

195 

piP+f^+xy'^) 

0.14 

0.09 

2.45 

1.10 

2.70 

2014 Jan 31 

TRAPPIST 

Occultation 

Sloan- 2 :' 

12 

947 

p(a^) + 0 

0.25 

0.10 

1.12 

1.00 

1.12 

2014 Feb 05 

TRAPPIST 

Occultation 

Sloan- 2 :' 

12 

1033 

p{P+xy^) + 0 

0.39 

0.17 

1.31 

1.49 

1.95 

2014 Mar 22 

TRAPPIST 

Occultation 

Sloan- 2 :' 

12 

1007 

p{P+xy^) + 0 

0.38 

0.16 

1.17 

1.18 

1.38 

2014 Apr 07 

TRAPPIST 

Transit 

Sloan- 2 :' 

13 

700 

+xy^) 

0.48 

0.19 

1.61 

1.45 

2.32 

2014 Apr 14 

TRAPPIST 

Occultation 

Sloan- 2 :' 

11 

851 

p(ai) 

0.36 

0.16 

1.02 

1.58 

1.61 

2014 Nov 08 

TRAPPIST 

Transit 

Johnson-S 

7 

966 

p(P+}p+xy~^) 

0.55 

0.23 

1.46 

1.09 

1.59 

2014 Dec 01 

EulerCam 

Transit 

Geneva-S 

90 

162 

p{a^+f^+xy^) 

0.11 

0.11 

1.51 

1.05 

1.59 

2014 Dec 24 

TRAPPIST 

Transit 

Sloan- 2 :' 

8 

961 

p(al) 

0.29 

0.12 

0.91 

2.07 

1.87 

2014 Dec 29 

EulerCam 

Transit 

Geneva-S 

60 

223 

P(a^+/^) 

0.08 

0.06 

1.15 

1.19 

1.36 


the average RV error, see Santerne et al. 2015) is 76.7 ms“^. 
These variations do not phase with the transit ephemeris, as 
one might expect if the observed RV signal was originating 
from a false-positive scenario, such as a blended eclipsing bi¬ 
nary (see e.g. Santos et al. 2002). However, as the scatter in 
the bisector span values is comparable to the semi-amplitude 
of the RV signal, we were not able to discard blend scenar¬ 
ios based on the traditional bisector span technique (Queloz 
et al. 2001). Instead, we performed a detailed blend analysis 
and used the FASTIS Bayesian software (Diaz et al. 2014, 
Santerne et al. 2014) to statistically validate the planet, as 
described in Section 4. 

2.3 Follow-up eclipse photometry 

To refine the system’s parameters, high-precision eclipse 
(transit and occultation) lightcurves were obtained using 
the 60cm TRAPPIST robotic telescope (TRAnsiting Plan¬ 
ets and Planetesimals Small Telescope) and the EulerCam 
CCD camera that is mounted on the 1.2m Euler-Swiss tele¬ 
scope, which are both located at ESO La Silla Observatory. 
These follow-up lightcurves are summarized in Table 2 and 
presented in Fig. 4 and 5. The transits were observed in dif¬ 
ferent filters to search for a potential color dependence of the 
transit depth, which might have been indicative of a blend 
(see Section 4). 

2.3.1 TRAPPIST observations 

TRAPPIST is a 60cm robotic telescope dedicated to the 
detection and characterization of transiting exoplanets and 
to the photometric monitoring of bright comets and other 
small bodies. It is equipped with a thermoelectrically-cooled 
2Kx2K CCD, which has a pixel scale of 0.65” that translates 
into a 22’x22’ field of view. For details of TRAPPIST, see 
Gillon et al. (2011) and Jehin et al. (2011). TRAPPIST was 


used to observe four transits of WASP-121 b in a Sloan- jz' 
filter (Aefi = 895 ± 1 nm) and one transit in a Johnson-R 
filter (Aefi = 440± 1 nm). As we noticed from a preliminary 
analysis that WASP-121 b is actually an extremely favorable 
target for secondary eclipse measurements (we will elaborate 
on this in Section 5.4), we also observed seven occultation 
windows in the Sloan- 2 ' filter. During the runs, the posi¬ 
tions of the stars on the chip were maintained to within a 
few pixels thanks to a “software guiding” system that reg¬ 
ularly derives an astrometric solution for the most recently 
acquired image and sends pointing corrections to the mount 
if needed. After a standard pre-reduction (bias, dark, and 
flatfield correction), the stellar fluxes were extracted from 
the images using the IRAF/DAOPHOT^ aperture photom¬ 
etry software (Stetson 1987). For each lightcurve, we tested 
several sets of reduction parameters and kept the one giving 
the most precise photometry for the stars of similar bright¬ 
ness as the target. After a careful selection of reference stars, 
the photometric lightcurves were finally obtained using dif¬ 
ferential photometry. 

2.3.2 EulerCam observations 

EulerCam is a 4Kx4K E2V CCD installed at the Cassegrain 
focus of the 1.2m Euler-Swiss telescope. The field of view of 
EulerCam is 15.7’x 15.7’, producing a pixel scale of 0.23”. 
To keep the stars on the same locations on the detector 
during the observations, EulerCam employs an “Absolute 
Tracking” system that is very similar to the one of TRAP¬ 
PIST, which matches the point sources in each image with 


^ IRAF is distributed by the National Optical Astronomy Obser¬ 
vatory, which is operated by the Association of Universities for 
Research in Astronomy, Inc., under cooperative agreement with 
the National Science Foundation. 
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Figure 4. Follow-up transit photometry for WASP-121 b. The 
observations are binned per 2 min and period-folded on the best- 
fit transit ephemeris (see Section 3.3). Each lightcurve has been 
divided by the respective best-fit photometric baseline model. For 
each filter, the superimposed, solid black line is our best-fit transit 
model. The lightcurves are shifted along the y-axis for clarity. 

a catalog, and if needed, adjusts the telescope pointing be¬ 
tween exposures to compensate for drifts. EulerCam was 
used to observe two transits of WASP-121 b in a Gunn-r' fil¬ 
ter (Aefi = 664±1 nm) and two other transits in a Geneva-i? 
filter^ (Aeff = 425 ± 1 nm). A slight defocus was applied 
to the telescope to optimize the observation efficiency and 
to minimize pixel-to-pixel effects. The reduction procedure 

^ http: //obswww. unige . ch/gcpd/phl3. html 
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Figure 5. To'p: Individual follow-up occultation lightcurves for 
WASP-121 b. The observations are binned per 2 min and period- 
folded on the best-fit transit ephemeris (see Section 3.3). Each 
lightcurve has been divided by the respective best-fit photomet¬ 
ric baseline model. For each lightcurve, the superimposed, solid 
black line is our best-fit occultation model. The lightcurves are 
shifted along the y-axis for clarity. Bottom: Combined follow-up 
occultation photometry for WASP-121 b (bin width=10 min). 
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Figure 6. Out-of-eclipse photometric monitoring: globally 
normalized TRAPPIST differential photometry obtained for 
WASP-121 in Johnson-_B (top), Johnson-V (middle), and Sloan- 
z' (bottom) filters, unbinned (cyan) and binned by day (black). 


used to extract the transit lightcurves was similar to that 
performed on TRAPPIST data. Farther details of the En- 
lerCam instrument and data reduction procedures can be 
found in Lendl et al. (2012). 


2.4 Out-of-eclipse photometric monitoring 

To search for potential out-of-eclipse photometric variability 
that would not have been detected in the WASP photome¬ 
try (see Section 2.1), we monitored WASP-121 with TRAP¬ 
PIST for 27 non-consecutive nights between 2014 Oct 25 and 
2014 Dec 8. This monitoring consisted in taking every night 
a short sequence of ten images in three filters: Johnson-R 
(Aeff = 440 ± 1 nm), Johnson-P (Aeff = 546.5 ± 1 nm), and 
Sloan-z' (Aeff = 895 ± 1 nm). The data were reduced as 
described in Section 2.3.1. The globally normalized differen¬ 
tial lightcurves obtained in each filter are shown in Fig. 6. 
WASP-121 appears to be very quiet in photometry, the stan¬ 
dard deviations of the binned lightcurves being 1.6 mmag 
(B), 1.3 mmag (V), and 1.1 mmag (z'). 


3 ANALYSIS 

3.1 Spectroscopic analysis 

The individual CORALIE spectra were co-added to produce 
a single spectrum with an average S/N of around 150:1. The 
analysis was performed using standard pipeline reduction 
products and the procedures given in Doyle et al. (2013). 
The derived stellar parameters are listed in Table 4. 

The excitation balance of the Fe l lines was used to de¬ 
termine the effective temperature Teff. The surface gravity 
log g* was determined from the ionisation balance of Fe l 
and Fe ll. The Ca l line at 6439 A and the Na l D lines were 
also used as log g* diagnostics. The iron abundance was 


determined from equivalent width measurements of several 
unblended lines and is relative to the solar value obtained by 
Asplund et al. (2009). A value for microturbulence (^t) was 
determined from Fe l using the method of Magain (1984). 
The quoted error estimates include those given by the un¬ 
certainties in Teff and log g*, as well as the scatter due to 
measurement and atomic data uncertainties. 

The projected stellar rotation velocity u* sin i* was de¬ 
termined by fitting the profiles of several unblended Fe l 
lines. A macroturbulent velocity (umac) of 6.0 ± 0.6 kms“^ 
was assumed using the asteroseismic-based calibration of 
Doyle et al. (2014) and an instrumental resolution of 55,000. 
A best-fitting value of u* sin i* = 13.5 ± 0.7 km s“^ was 
obtained. 

There is no significant detection of lithium in the spec¬ 
tra, with an abundance upper limit log A(Li)<1.0. The lack 
of any detectable lithium does not provide an age constraint 
as the star’s T^s places it in the lithium gap (Bohm-Vitense 
2004). There is also no significant chromospheric emission 
in the Ca ll H and K line cores. 

The spectral type was estimated from Teff using the 
Table B.l in Gray (2008) and the Torres et al. (2010) cali¬ 
bration was used to obtain first stellar mass and radius es¬ 
timates: M* = 1.37 ± 0.14 Mq and R* = 1.52 ± 0.41 Rq. 


3.2 Stellar jitter 

As mentioned in Section 2.2, WASP-121 exhibits an espe¬ 
cially high scatter in its RV residuals (Fig. 2a), CCF bisector 
spans (Fig. 2c), and FWHM (Fig. 2d). Fig. 7 compares the 
RV residuals to the CCF bisector spans assuming a circular 
(top) and an eccentric (bottom) orbit. There is a significant 
anti-correlation between these two quantities, the correla¬ 
tion coefficients being -0.68 and -0.67 in the circular and 
eccentric cases, respectively. Such an anti-correlation is com¬ 
monly interpreted as being a signature of stellar activity (see 
e.g. Queloz et al. 2001, Melo et al. 2007), but Santerne et al. 
(2015) showed that it could also be produced by a blended 
star with a lower CCF FWHM (and thus slower rotation) 
than the target star. However, as detailed in Section 4, we 
were not able to reproduce the observed RVs, CCF bisector 
spans, and FWHM assuming such a scenario, thus making 
it likely that the high jitter measured for WASP-121 is due 
to stellar activity. One might have concerns about the non¬ 
detection of this stellar activity in the form of emission in 
the Ca II H and K line cores (see Section 3.1), out-of-eclipse 
photometric variability (see Sections 2.1 and 2.4), or spot¬ 
crossing events during transits (see Fig. 4). A detailed study 
of WASP-121’s activity being beyond the scope of this pa¬ 
per, we just propose here some potential explanations re¬ 
garding these non-detections. 

First, we note that such a situation is not atypical 
for an F-type star as it was also encountered for the exo¬ 
planet F-type host stars HAT-P-33 (Hartman et al. 2011) 
and WASP-111 (Anderson et al. 2014), which both present 
high activity-related RV jitter with no other apparent sign 
of stellar activity. More precisely, in the case of WASP-111, 
the scatter in the RV residuals, CCF bisector spans and 
FWHM dropped from one season to the next, identifying 
clearly stellar activity as an origin for these. 

Secondly, it can be seen from Noyes et al. (1984) that 
the chromospheric Ca ll emission of stars decreases with low- 
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CCK BIS [m/s 



CCF BIS [rri/s] 


Figure 7. RV residuals from the best-fit circular (top) and ec¬ 
centric (bottom) orbital models vs. CCF bisector spans. The two 
quantities are clearly anti-correlated. Linear fits to the data are 
overplotted and the correlation coefficients (r) are given. 


ering B — V. As WASP-121’s B — V is only 0.43, its Ca ii 
emission could simply be too weak to be detected in our low 
S/N CORALIE spectra. 

As for the non-detections of out-of-eclipse variations 
or spot-crossing events in the photometry, they might be 
explained, at least to some extent, if the star is plage- 
dominated, in opposition to spot-dominated (see Dumusque 
et al. 2014). Both spots and plages are regions of strong 
magnetic fields that inhibit locally the convection and sup¬ 
press the convective blueshift effect (see e.g. Dravins et al. 

1981) . These regions are thus redshifted compared to the 
quiet photosphere and induce RV variations as they appear 
and disappear from the visible stellar disk due to the rota¬ 
tion of the star. The amplitude of this effect is expected to 
be similar for spots and plages of the same size. However, 
plages present a much lower flux ratio with the quiet pho¬ 
tosphere compared to spots (see Dumusque et al. 2014 for 
details), so that they induce smaller brightness variations as 
they appear and disappear from the visible stellar disk. 

Fig. 8 shows Lomb-Scargle periodograms (Scargle 

1982) of the RV residuals (assuming a circular orbit, second 
panel from the top), CCF bisector spans (third panel from 
the top), and FWHM (bottom panel). In each of these 
periodograms, the highest peaks are found at periods of 
0.89, 1.13, and 8.3 days (all with false alarm probabilities 
< 0.001, except for the peak at 1.13 days in the periodogram 
of the bisector spans). We could assume that one of these 
three periods, which are each other’s daily aliases, is the 
stellar rotation period (this supposes that the activity- 
induced RV signal is related to the rotation of the star). A 
stellar rotation period of 8.3 days does not seem very likely, 
as the maximal rotation period implied by w* sin i* and our 
final estimate of the stellar radius (R* = 1.458 ± 0.030 R©, 
see Section 3.3) is 5.46 ± 0.32 days (assuming sin i* = 1). 
Stellar rotation periods of 0.89 and 1.13 days would imply, 
together with our measured value for R*, rotation velocities 



1 10 100 
Period [d] 




Oh )_ 



Period [d] 


Figure 8. From top to bottom: Lomb-Scargle periodograms of 
the RVs, RV residuals (assuming a circular orbit), CCF bisector 
spans, and FWHM. Horizontal lines indicate false alarm proba¬ 
bility levels of 0.001, 0.01, and 0.1 (from top to bottom). 


Vi, of 82.9 ± 1.7 kms“^ and 65.3 ± 1.4 kms“^, respectively, 
and inclinations i* of 9.4 ± 0.6 deg and 11.9 ± 0.7 deg, 
respectively, where the uncertainties are the quadratic 
sum of the uncertainties due to each input parameter. In 
both cases, the star would thus be seen nearly pole-on, i.e. 
with the rotation axis oriented towards the Earth. This 
configuration could also contribute to the non-detection of 
out-of-eclipse photometric variations, as the modulation 
in brightness due to the changing visibility of star spots 
or plages with rotation could then be negligible. If the 
star is seen nearly pole-on, WASP-121b should then be in 
a (nearly) polar orbit to produce transits. As mentioned 
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previously, this is actually what we hnd via our observations 
of the RM effect (see Sections 3.3 and 5.2), which makes the 
stellar rotation period much likely to be 0.89 or 1.13 days. 
We note that for the 1-Gyr-old cluster NGC6811 observed 
with Kepler by Meibom et al. (2013), the stars with 
B — V = 0.43 (such as WASP-121) have ~l-day rotation 
periods (see their Fig. lb). As the age of WASP-121 
obtained by stellar evolution modeling is 1.5 ± 1.0 Gyr 
(see Section 3.4), a ~l-day rotation period would thus be 
consistent with their values. 

3.3 Global modeling of the data 

To determine the system parameters, we performed a com¬ 
bined analysis of the follow-up photometry and the RV data, 
using for this purpose the adaptive Markov Chain Monte- 
Carlo (MCMC) code described in Gillon et al. (2012) and 
references therein. To model the photometry, we used the 
eclipse model of Mandel & Agol (2002) multiplied by a differ¬ 
ent baseline model for each lightcurve. These baseline mod¬ 
els aim to represent astrophysical, instrumental, or environ¬ 
mental effects, which are able to produce photometric varia¬ 
tions and can, therefore, affect the photometric lightcurves. 
They can be made up of different first to fourth-order poly¬ 
nomials with respect to time or other variables, such as air- 
mass, PSF full-width at half maximum, background, or stel¬ 
lar position on the detector. To find the optimal baseline 
function for each lightcurve, i.e. the model minimizing the 
number of parameters and the level of noise in the best-fit 
residuals, the Bayes factor, as estimated from the Bayesian 
Information Criterion (BIG, Schwarz 1978), was used. The 
best photometric baseline functions are listed in Table 2. For 
eight TRAPPIST lightcurves (see Table 2), a normalization 
offset was also part of the baseline model to represent the 
effect of the meridian flip; that is, the 180° rotation that the 
German equatorial mount telescope has to undergo when 
the meridian is reached. This movement results in different 
positions of the stellar images on the detector before and 
after the flip, and the normalization offset allows to take 
into account a possible consecutive jump in the differential 
photometry at the time of the flip. 

On their side, the RVs were modelled using a classi¬ 
cal Keplerian model (e.g. Murray & Correia 2010) added to 
a baseline model for the stellar and instrumental variabil¬ 
ity. For the RVs obtained during a transit, the RM effect 
was modelled using the formulation of Gimenez (2006). The 
RVs were partitioned into four datasets, each with a differ¬ 
ent baseline model: the RVs obtained before the replacement 
of the CORALIE optical fibre (37 RVs, dataset #1), those 
obtained after (15 RVs, dataset #2), the first RM sequence 
(19 RVs, dataset #3), and the second one (18 RVs, dataset 
#4). For all datasets, the minimal baseline model was a 
scalar V-y representing the systemic velocity of the star. For 
datasets and ^2, first-order polynomial functions of the 
CCF bisector span and FWHM were also part of the baseline 
models. This choice of baselines allowed to reduce the scat¬ 
ter in the global RV residuals from 67.1 ms“^ to 37.7 ms“^ 
(see Fig. 2a and b) and was thus strongly favored by the 
BIG. These additional baseline terms were not necessary for 
the datasets #3 and ^4 (RM sequences). 

The basic jump parameters in our MCMC analy¬ 
ses, i.e. the parameters that are randomly perturbed 


Filter 

Ui 

U2 

Sloan-^' 

0.171 ± 0.014 

0.306 ± 0.006 

Gunn-r' 

0.295 ± 0.015 

0.325 ± 0.006 

Johnson- and Geneva-B 

0.510 ± 0.027 

0.260 ± 0.019 


Table 3. Expectations and standard deviations of the normal 
distributions used as prior distributions for the quadratic limb- 
darkening (LD) coefficients ui and U 2 in our MCMC analysis. 

at each step of the MCMC, were: the planet/star area 
ratio dF = (Rp/R*)^; the occultation depth in the 2 :'-band 
dRocc, z '; the transit impact parameter in the case of a circu¬ 
lar orbit b' = a cos ip/Ry,, where a is the orbital semi-major 
axis and ip is the orbital inclination; the transit width (from 
1st to 4th contact) W ; the time of mid-transit Tq; the orbital 
period R; the stellar effective temperature Teff; the stellar 
metallicity [Fe/H]; the parameter K 2 = Ky/1 — R^^®, 

where K is the RV orbital semi-amplitude and e is the 
orbital eccentricity; the two parameters y/e cos oj and 
^/e sin uj, where uj is the argument of the periastron; and 
the two parameters y/vy, sin i* cos j3 and V 'v* sin i* sin j3, 
where u* sin i* is the projected rotational velocity of the 
star and /? is the sky-projected angle between the stellar 
spin axis and the planet’s orbital axis. The reasons to use 
y/e cos OJ and y/e sin ut as jump parameters instead of the 
more traditional e cos oj and e sin oj are detailed in Triaud 
et al. (2011). The results obtained from the spectroscopic 
analysis (see Section 3.1) were used to impose normal 
prior distributions on Teff, [Fe/H], and u* sin with 
expectations and standard deviations corresponding to the 
quoted measurements and errors, respectively. Uniform 
prior distributions were assumed for the other parameters. 
The photometric and RV baseline model parameters were 
not actual jump parameters; they were determined by 
least-square minimization from the residuals at each step of 
the MCMC. 

The effect of stellar limb-darkening on our transit 
lightcurves was accounted for using a quadratic limb- 
darkening law, where the quadratic coefficients ui and U 2 
were allowed to float in our MCMC analysis. However, 
we did not use these coefficients themselves as jump 
parameters but their combinations, Ci = 2 x Ui -|- M 2 and 
C 2 = Ml — 2 X U 2 , to minimize the correlation of the obtained 
uncertainties as introduced by Holman et al. (2006). To 
obtain a limb-darkening solution consistent with theory, 
we used normal prior distributions for mi and M 2 based on 
theoretical values and l-cr errors interpolated in the tables 
by Claret & Bloemen (2011). These prior distributions are 
presented in Table 3. 

A preliminary analysis was performed to determine 
the correction factors {CF) for our photometric errors, as 
described in Gillon et al. (2012). For each lightcurve, CF 
is the product of two contributions, P-w and Pr- On one 
side, P-w represents the under- or overestimation of the 
white noise of each measurement. It is computed as the 
ratio between the standard deviation of the residuals and 
the mean photometric error. On the other side, Pr allows 
us to take into account the correlated noise present in the 
lightcurve (i.e., the inability of our model to perfectly fit 
the data). It is calculated from the standard deviations of 
the binned and unbinned residuals for different binning 
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intervals ranging from 5 to 120 min with the largest value 
being kept as fir- The standard deviation of the best-ht 
residuals (unbinned and binned per intervals of 2 min) 
and the deduced values for /J™, /3r and CF = l3w x /3r for 
each lightcurve are presented in Table 2. This preliminary 
analysis also allowed to compute the jitter values that were 
added quadratically to the RV errors of each RV dataset 
to equal their mean values to the standard deviations of 
the best-fits residuals, and thus achieve reduced values 
of unity. These jitter values were 27.9 ms“^, 22.3 ms“^, 
16.9 ms“^, and zero for the datasets #2, #3, and ^4, 
respectively. 

A second analysis with the updated photometric and 
RV errors was then performed to determine the stellar 
density p*, which can be derived from the Kepler’s third 
law and the jump parameters (Rp/R*)^, b', W, P, y/e cos tu 
and ^/e sin u) (see e.g. Seager & Mallen-Ornelas 2003 and 
Winn 2010). This analysis consisted of two Markov chains 
of 10® steps, whose convergence was checked using the 
statistical test of Gelman & Rubin (1992). The hrst 20% 
of each chain was considered as its burn-in phase and 
discarded. The resulting stellar density was used as input 
of a stellar evolution modeling, together with the effective 
temperature and metallicity derived from spectroscopy, 
with the aim to estimate the stellar mass M* and the age 
of the system. This stellar evolution modeling is described 
in details in Section 3.4. 

Two final analyses were then performed: one assuming 
a circular orbit (e = 0) and one with a free eccentricity. 
Each analysis consisted of two Markov chains of 10® steps, 
whose convergence was again checked with the Gelman & 
Rubin test (Gelman & Rubin 1992). As previously, the 
first 20% of each chain was considered as its burn-in phase 
and discarded. At each step of the Markov chains, p* was 
computed as described above and a value for M* was 
drawn within a normal distribution having as expectation 
and standard deviation the value and error obtained from 
the stellar evolution modeling. The stellar radius R* was 
derived from M* and p*, and the other physical parameters 
of the system were then deduced from the jump parameters 
and stellar mass and radius. It appeared that the circular 
orbit was strongly favored by the Bayes factor (~5000 in 
its favor) compared to the eccentric orbit. As there was no 
evidence for a significant eccentricity (e < 0.07 at 3 ct), we 
thus adopted the circular orbit as our nominal solution. The 
corresponding derived parameters are presented in Table 4, 
while the best-fit models are shown in Fig. 2b (RVs), 3 
(zoom on the RM effect), 4 (transit photometry), and 5 
(occultation photometry). 

3.4 Stellar evolution modeling 

As introduced in Section 3.3, we performed for the host 
star a stellar evolution modeling based on the CLES code 
(Scuflaire et al. 2008), in order to estimate the stellar mass 
and the age of the system. We used as inputs the stel¬ 
lar density deduced from the transit lightcurves, and the 
effective temperature and metallicity derived from spec¬ 
troscopy (see Table 4). We considered that [Fe/H] repre¬ 
sents the global metallicity with respect to the Sun i.e. 
[Fe/H] = [log(%/X)* - log(Z/X)0], with {Z/X)q = 0.0181 
(Asplund et al. 2009). The parameter of the mixing-length 


WASP-121 



TeflfKl 

Figure 9. Evolutionary tracks in a — pi,lpQ HR diagram for 
WASP-121, for different masses and metallicities, within (solid 
lines) and without (dashed lines) the l-cr box Tefj — p*/p 0 . Some 
stellar ages are also indicated. 


theory (MLT) of convection was kept fixed to the so¬ 
lar calibration (omct = 1.8), and the possible convective 
core extra-mixing (due to overshooting, rotationally-induced 
mixing, etc.) and microscopic diffusion (gravitational set¬ 
tling) of elements were included. 

We obtained a stellar mass of 1.355 ± 0.080 Mq and 
an age of 1.5 ± 1.0 Gyr. These l-cr uncertainties were ob¬ 
tained by considering the respective l-cr range for the effec¬ 
tive temperature, metallicity and stellar density, but also by 
varying the internal stellar physics. We computed, since the 
helium atmospheric abundance cannot be directly measured 
from spectroscopy, evolutionary tracks with two initial he¬ 
lium abundances: the solar value (F/i.o = 0.27), and a value 
labelled Yq that increases with Z (as expected if the local 
medium follows the general trend observed for the chemical 
evolution of galaxies; Izotov & Thuan 2010). We also in¬ 
vestigated the effects of the possible convective core extra¬ 
mixing, by varying it between 0 and 0.3 (Noels & Montalban 
2013). 

Evolutionary tracks are presented on Fig. 9 for sev¬ 
eral stellar masses and metallicities. Obtaining an accurate 
stellar mass from evolution modeling primarily needs accu¬ 
rate spectroscopic estimates for the effective temperature 
but also, very importantly, for the metallicity. Metallicity 
is a parameter that strongly determines the location of the 
evolutionary tracks in a HR diagram (compare in Fig. 9 the 
solid and dashed tracks for identical stellar masses, but with 
different metallicities). The stellar mass obtained from stel¬ 
lar evolution modeling (1.355 ± 0.085 Mq) is in excellent 
agreement with the first estimate derived from the Torres 
et al. (2010) calibration in Section 3.1. The inferred stel¬ 
lar age of 1.5 ± 1.0 Gyr places WASP-121 on the main se¬ 
quence. The total lifetime on the main sequence (H-core 
burning) of WASP-121 is 4.3 Gyr for a moderate extra¬ 
mixing aem = 0.2, and is reduced to 3.3 Gyr without con¬ 
sidering any extra-mixing process in the core. 


MNRAS 000, 1-19 (2014) 





10 L. Delrez et al. 


General information 

RA (J2000) 

07h 10m 24.06s 

U-magnitude 

10.44 

Dec (J2000) 

-39° 05' 50.55" 

A-magnitude 

9.37 

Distance [pc]“ 

270 ± 90 



Stellar parameters from spectroscopic analysis 

Spectral type 

F6V 

Microturbulence [km s~^] 

1.5 ± 0.1 

Effective temperature Teff [K] 

6460 ± 140 

Macroturbulence Umac [km s“^] 

6.6 ± 0.6 

Surface gravity log [cgs] 

4.2 ± 0.2 

Proj. rot. velocity u* sin [km s“^] 

13.5 ± 0.7 

Metallicity [Fe/H] [dex] 

0.13 ± 0.09 

Lithium abundance log A(Li) [dex] 

< 1.0 


Parameters from MCMC analysis 


Jump parameters 


Planet/star area ratio (i?p/i?*)^ [%] 

1.551 ± 0.012 

Effective temperature Teff [K]^ 

6459 ± 140 

Occultation depth [ppm] 

603 ± 130 

Metallicity [Fe/H] [dex]® 

0.13 ± 0.09 

&'= a cos 2p/i?* [R-k] 

U.iDU_Q Q42 

Cl,2' = 2 111,2/ + “2,z' 

0 n‘:>7-HU.uzb 

U.DO/_o.o25 

Transit width W [d] 

0.1203 ± 0.0003 

C2,z' = “l,z' - 2 ll2,z' 

-0.445 ± 0.018 

To - 2 450 000 [HJDtdb] 

6635.70832l°;00““ 

Cj^ y! — 2 12f -j- 222 t' 

u.yu^_Q Q26 

Orbital period P [d] 

1 974Q9c;c:+0-00000020 

QQQQQQ25 

— ^l,r^ ^ ^2,r^ 

—0 351+0-020 
U.OOi_Q Qlg 

RV K2 [ms-i di/3] 

196.4l®:« 

ci,fl = 2 111,a + 112,B 

i.zuu_Q 63;,^ 

y/e cos iJ 

0 (fixed) 

C2,S = iii,s — 2 112,s 

u.uiz_o 045 

y/e sin (jj 

0 (fixed) 



y/v-k sin 2* cos (3 

-0.78 ±0.34 



y/Vi, sin 2* sin {3 

—3 



Deduced stellar parameters 

Mean density p* [pq] 

n 497+0008 
U.4c5/_o.o09 

Limb-darkening coefficient Ui 

0.166 ± 0.013 

Surface gravity log p* [cgs] 

4.Z4Z_o_oi2 

Limb-darkening coefficient U 2 z' 

0.305 ± 0.008 

Mass A4* [M©]° 

1 ocq-I-0.080 

1.009_q 079 

Limb-darkening coefficient Ui ^1 

0.290 ± 0.014 

Radius R* [Rq] 

1.458 ± 0.030 

Limb-darkening coefficient U 2 r' 

0.325 ± 0.007 

Luminosity L* [-^^0] 

3.3 ± 0.3 

Limb-darkening coefficient ui b 

0.505 ± 0.018 

Proj. rot. velocity u* sin 2* [km s~^]^ 

10.ou_o.68 

Limb-darkening coefficient 1x2,s 

U.ZOU_Q Qj^g 

Mean systemic velocity Vy [kms~^] 

38.350 ± 0.021 



Deduced planet parameters 

RV K [ms-i] 

181.ll®i 

Mean density pp [pjup] 

0.201 ± 0.010 

Planet/star radius ratio Rp/R* 

0.12454l°;00°^^ 

Surface gravity log pp [cgs] 

2.973 ± 0.017 

Tocc - 2 450 000 [HJDtdb] 

6636.345781°;“““ 

Mass Mp [Mjup] 

1 1 oq+0-064 

i.i09_Q Q62 

Scaled semi-major axis a/R^ 

0 7C4-I-O.O23 

9. /O4_0.o28 

Radius Rp [Rjup] 

1.807 ± 0.039 

Orbital semi-major axis a [AU] 

'^•^^^^^-0.00050 

Roche limit ur [AU]^ 

0.02205 ± 0.00066 

Orbital inclination 2p [deg] 

87.6 ± 0.6 

a/uR 

1.153 ± 0.019 

Orbital eccentricity e 

0 (fixed) 

Equilibrium temperature Teq [K]® 

2358 ± 52 

Argument of periastron oj [deg] 

- 

Irradiation [erg s“^cm~^] 

7.lll;° 10° 

Sky-projected orbital obliquity I3 [deg] 

257.8l^;| 



Planet parameters corrected for asphericity 

Radius Rp [Rjup] 

1.865 ± 0.044 

Mean density pp [pjup] 

0.183 ± 0.016 


Table 4. System parameters for WASP-121. The values given for the parameters derived from our MCMC analysis are medians and 
l-(7 limits of the marginalized posterior probability distributions. “From V mag and estimated absolute magnitude. ^Using as pri¬ 
ors the values derived from the spectroscopic analysis. “Using as prior the value obtained from stellar evolution modeling. ^Using 
Ur = 2.46 (Chandrasekhar 1987). “Assuming a null Bond albedo and isotropic reradiation (reradiation factor /=l/4, 

Lopez-Morales &; Seager 2007). 
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4 FASTIS VALIDATION 

Astrophysical false-positive scenarios such as eclipsing bina¬ 
ries might mimic both the transit and radial velocity signal 
of a planet (Torres et al. 2005). In some particular configu¬ 
rations, even the Rossiter-McLaughlin effect might be mim¬ 
icked (Santerne et al. 2015). In the present case, the line 
bisector and the FWHM present a large variation. Moreover 
we were able to detect from the ground in the a'-band a rel¬ 
atively deep secondary eclipse. Both arguments prevent to 
secure the planetary nature of WASP-121 b without a more 
careful investigation. We can imagine five scenarios which 
can reasonably produce the observed data: (1) a planet tran¬ 
siting the target star, (2) a planet transiting a chance-aligned 
background star, (3) a chance-aligned background eclipsing 
binary, (4) a planet transiting a physical companion to the 
target star, and (5) a star eclipsing a physical companion 
to the target star (i.e. a triple system). The scenarios 2-5 
can be split in two categories and discussed separately: the 
background or physical companion sources of false positive. 

4.1 A background source 

A background source of false positive can mimic the transit 
data of a planet within a range a magnitude Am defined as 
(Morton & Johnson 2011): 

Am = 2.51ogio , (1) 

where 5tr and Sf,g are the depth of the transit, as measured 
in the lightcurve, and the depth of the background eclipse, 
respectively. Assuming a maximum eclipse depth of 50% for 
the background star, we find that the maximum magnitude 
range is of 3.78. Note that a 50% depth eclipsing binary 
should also produce a secondary eclipse with a similar depth, 
which is clearly excluded in the case of WASP-121. There¬ 
fore, this magnitude range is over-estimated, which will also 
over-estimate the probability of a background star as the 
host of the observed signal. The star is of magnitude 10.44 
in the V-band, hence false positive can probe stars as faint 
as magnitude 14.22 in the same bandpass. 

To compute the probability of having a background star 
chance-aligned with WASP-121, we took the APASS (Hen- 
den et al. 2015) DR8 V-magnitude of all stars within 1 de¬ 
gree from WASP-121. The magnitude limit of this catalog 
Is between magnitude 16 and 17, which is about 2 mag¬ 
nitudes above the maximum magnitude of the background 
star. We therefore assume that this catalog Is complete in 
the range of magnitude we are considering here. We did not 
detect in the TRAPPIST data any background source that 
could host the transit. Thus, in this scenario, the background 
source should be blended within the TRAPPIST PSF. In 
Fig. 10 we display the sensitivity curve of the TRAPPIST 
PSF. This curve was obtained by injecting artificial stars in 
good-seeing TRAPPIST images at various separations and 
with a range of magnitude, and then attempting to detect 
them with IRAF/DAOPHOT. Any star brighter than mag¬ 
nitude 14 and separated from WASP-121 by more than 3.5” 
should have been detected. Using these constraints and as¬ 
suming that the stars are randomly distributed around the 
target star, we find that the probability for WASP-121 to be 
chance-aligned with a background star brighter than magni¬ 
tude 14.22 is at the level of maximum 8.10“^. If we account 





Heliocentric radial velocity [km.s^^] 


Figure 10. Upper panel: Map of the density of background stars 
chance-aligned with WASP-121, integrated within an angular sep¬ 
aration of up to 10”, as function of the U-band magnitude. The 
positive-slope hatched region displays all the stars that would 
have been significantly detected in the TRAPPIST data. The 
negative-slope hatched region displays the stars that are too faint 
to reproduce the observed transit depth of WASP-121. Lower 
panel: Distribution of heliocentric radial velocity of stars within 
10 degrees from WASP-121 as observed by the RAVE experi¬ 
ment. The hatched regions show the velocity of background stars 
that should be spectroscopically resolved with WASP-121 and 
thus should produce a radial velocity variation in anti-phase with 
the transit ephemeris, or no significant variation (Santerne et al. 
2015). 


for the probability that this hypothetical background star 
host a binary (~50% ; Raghavan et al. 2010) or a planet 
(~50% ; Mayor et al. 2011) and the eclipse or transit proba¬ 
bility at 1.3 days (~25%), we end with a a priori probability 
that the transit signal is hosted by a background star lower 
than ~ 1.10“"^. 

A significant radial velocity variation has been detected 
in phase with the transit ephemeris and no other component 
has been significantly detected in the CORALIE CCFs with 
a flux ratio greater than 0.7% (i.e. 5.4 magnitudes). There¬ 
fore, to mimic both the photometric and radial velocity data 
of WASP-121, the background source should be chance- 
aligned along the line of sight and in the radial velocity 
space. To produce a radial velocity in phase with the transit 
ephemeris, the background star should have a systemic ve¬ 
locity within the width of the line profile of WASP-121 (San¬ 
terne et al. 2015), i.e. the two stars should be spectroscopi¬ 
cally unresolved. Given the systemic velocity of WASP-121 
(7 ~ 38.3 km s“^) and the width of its line profile observed 
by CORALIE {accF = FWHM/2^2ln(2) « 8.2 km s“^), 
the background star should thus have a systemic radial ve¬ 
locity in the range [30.1 ; 46.5] km s“^. We then compared 
this range of radial velocity with the distribution of stars in 
the vicinity of WASP-121 from the RAVE database (Kor- 
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dopatis et al. 2013). The offset between the CORALIE zero- 
point and the RAVE reference is expected to be less than a 
few hundreds of m s“^. We took all the stars observed by 
this spectroscopic survey within 10 degrees of WASP-121 for 
which we display their radial velocity distribution in Fig. 10. 
Only 20% of the stars in the vicinity of WASP-121 should 
be spectroscopically unresolved. Assuming that the position 
of a star in a place of the sky is independent from its radial 
velocity, we find that the a priori probability that the signal 
observed in WASP-121 is caused by a background source is, 
at most, at the level of 20 ppm. We can therefore exclude 
all background false-positive scenarios as the source of the 
detected signal. 


4.2 A physical companion 

Apart from the planet scenario, the most likely scenario to 
reproduce both the photometric and spectroscopic data is a 
system physically bound with the target star. In such con¬ 
figuration, the various stellar components could be easily 
blended in both the plane of the sky and the radial veloc¬ 
ity space. Then, the transit could be mimicked either by an 
eclipsing star or a transiting sub-stellar object. 

To estimate the probability that WASP-121 is either a 
planet transiting the target, a physical companion (planet 
in binary) or a triple system, we used the FASTIS soft¬ 
ware (Diaz et al. 2014, Santerne et al. 2014, Santerne et al. 
2015) to model the llk+ photometric measurements ob¬ 
tained by TRAPPIST and EulerCam in different filters. 
The lightcurve was modelled using the EBOP code (Nelson 
& Davis 1972, Etzel 1981, Popper & Etzel 1981) extracted 
from the JKTEBOP package (Southworth 2008). For the limb 
darkening coefficients, we used the interpolated values from 
Claret & Bloemen (2011). To model the stars, we used the 
Dartmouth stellar evolution tracks of Dotter et al. (2008) 
and the BT-SETTL stellar atmosphere models of Allard 
et al. (2012) that we integrated in each individual band¬ 
pass. Since they are gravitationally bound, all stars were 
assumed to have the same metallicity and the same age. 
The orbits were assumed to be circular. We imposed that 
the physical companion is fainter by at least 1 magnitude in 
the V-band than the target star, otherwise, it would have 
been clearly identified in the spectral or CCF analysis. Each 
of the 16 lightcurves was modelled self-consistently allowing 
a different out-of-transit flux, contamination and an extra 
source of white noise (jitter). As FASTIS is not yet able to 
model the activity of stars in radial velocity data, we de¬ 
cided not to use these extra constraints. We also modelled 
the spectral energy distribution of WASP-121 composed by 
the magnitudes in the Johnson-B and -V, Sloan-^', -r', and 
-i', 2-MASS J, H, and As, and WISE W1 to W4 band- 
passes from the APASS database (Henden et al. 2015) and 
the AllWISE catalog (Wright et al. 2010). 

We analysed the aforementioned data using a MCMC 
procedure described in Dfaz et al. (2014). For the priors, we 
used the results from the spectroscopic analysis (Section 3.1) 
for the parameters of the target star and the initial mass 
function from Kroupa (2001) for the blended stars. For 
the orbital ephemeris, we used Normal priors matching the 
ephemeris reported in Table 4 with uncertainties boosted 
by 100, to avoid biasing the results with too narrow priors. 
For the other parameters, we chose uninformative priors. 
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Figure 11. Phase-folded binned photometric data of WASP-121 
observed in the various bandpasses together with the best planet 
(solid line), triple system (dash line), and planet in binary (dot 
line) models. 


We limited the priors on the planet radius to be less than 
2.2 i?jup, which is the radius of the biggest planet found so 
far: KOT13 (Szabo et al. 2011). The exhaustive list of pa¬ 
rameters and their priors are reported in Table Al. Both 
the planet and triple scenarios were described by ten free 
parameters, while the planet in binary scenario used 11 free 
parameters. Among them, eight were in common (target 
star and orbital parameters). An additional 49 free parame¬ 
ters were needed to describe the data: three for each of the 
16 lightcurves as already mentioned and an extra-source of 
white noise for the spectral energy distribution. For all sce¬ 
narios, we ran 40 chains of 3.10® iterations, randomly started 
from the joint prior distribution. All chains converged to¬ 
ward the same maximum of likelihood. We then thinned and 
merged the chains to derive the posterior distributions of all 
scenarios. They all have a posterior distribution with more 
than 1000 independent samples. We report in Table A2 the 
median and 68.3% confidence interval for the free param¬ 
eters. All the fitted parameters for the planet scenario are 
compatible within 1-a with those derived in Section 3.3. 

We display in Fig. 11 the transit and occultation data 
of WASP-121 in the various bandpasses together with the 
best planet, planet in binary and triple models. While the 
planet scenario reproduces well all the data, the triple sys¬ 
tem scenario is not able to fit perfectly the relatively sharp 
ingress and egress of the transit. The main difference be¬ 
tween the two scenarios is for the secondary eclipse. While 
the planet scenario is able to reproduce well the occultation 
data, the best triple scenario exhibits a secondary depth of 
more than 2000 ppm, which is clearly excluded by the data. 
For the planet transiting a companion star to the target, 
it is not possible to reproduce the observed transit depth 
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Figure 12. Left: Orbital distance to the Roche limit ratio-mass diagram for the known transiting hot Jupiters with 
0.2 Mjup < Mp < 13 Mjup and P < 12 d (data from the NASA Exoplanet Archive). We only show planets with a/a^ < 3. The 
planet WASP-121 b is shown in red. Right: Irradiation-radius diagram for the same sample of planets. 


and duration if we assume both that the companion star is 
fainter by at least 1 magnitude than the target star and that 
the maximum planetary radius allowed is 2.2 Pjup. In these 
conditions, to reproduce well the transit data, the planet 
would need to have a radius of 3.24 i?jup, which is clearly 
unphysical for such objects. 

We estimated the Bayes factor, Bij, between the three 
scenarios using the method of Tuomi & Jones (2012). This 
method has some limitations that are discussed in Santerne 
et al. (2014). However, since we tested here scenarios that 
have nearly the same number of free parameters and most 
of these parameters have the same priors in the various sce¬ 
narios, we assume that these limitations have no significant 
impact on our results. We found that logj^p Bij ~ 250 and 
logj^g Bij ~ 870, against the triple and planet-in-binary sce¬ 
narios (respectively), in favour to the planet one. This signif¬ 
icantly validates WASP-121 b as a bona-fide planet. Com¬ 
puting the BIC between the best models give a value of 
1103 and 3905 against the triple and planet-in-binary sce¬ 
narios (respectively), which confirms that WASP-121 b is a 
transiting planet. 

5 DISCUSSION 

WASP-121 b is a ~1.18 Mjup planet in a 1.27 days or¬ 
bit around a bright (U=10.4) F6V star. Its most notable 
property is that its orbital semi-major axis is only ~1.15 
times larger than its Roche limit, which suggests that the 
planet might be close to tidal disruption. Fig. 12 (left) shows 
the distribution of the orbital distance to the Roche limit 
ratio (a/au) as a function of the planetary mass for the 
known transiting planets with 0.2 Mjup < Mp <13 Mjup 
and P < 12 d. WASP-121 b is one of the closest sys¬ 
tems to tidal disruption, its direct competitors being 
WASP-12b (Hebb et al. 2009), WASP-19b (Hebb et al. 


2010), WASP-103b (Gillon et al. 2014), OGLE-TR-56b 
(Konacki et al. 2003, Adams et al. 2011), and WTS-2b 
(Birkby et al. 2014). According to Matsumura et al. (2010), 
these extreme planets are now expected to undergo tidal or¬ 
bital decay through tidal dissipation inside the star only 
(we will elaborate on this in Section 5.2 in the case of 
WASP-121 b). The speed of this final orbital decay depends 
mainly on the mass of the planet and the tidal dissipation 
efficiency of the star. As noted by Gillon et al. (2014), the 
fact that all the planets near tidal disruption found to date 
have similar masses (between 1.1 and 1.5 Mjup, see the left 
panel of Fig. 12) could thus suggest a narrow range of tidal 
dissipation efficiencies for the host stars of the known tran¬ 
siting hot Jupiters (mainly solar-type stars). 

5.1 Structural parameters of WASP-121 b: 
correction for asphericity using Roche 
geometry 

Being close to its Roche limit, WASP-121 b might be 
signihcantly deformed by the intense tidal forces it is 
subject to (e.g. Budaj 2011) and even lose some of its 
mass through Roche lobe overflow (e.g. Li et al. 2010, Lai 
et al. 2010). To evaluate the tidal distortion of the planet, 
we calculate its Roche shape by using the Roche model 
of Budaj (2011), as done by Southworth et al. (2015) for 
WASP-103 b. This model assumes that the planet is on a 
circular orbit and that it is rotating synchronously with 
its orbital period (Prot = Porb). This second assumption is 
perfectly reasonable as current theories of tidal evolution 
of close-in exoplanets predict synchronization times much 
shorter than times needed for circularization of the orbits 
(see e.g. Rasio & Ford 1996). The model takes as main 
inputs the orbital semi-major axis (a = 5.47 ± 0.11 P©), 
the star/planet mass ratio (M*/Mp = llOSll^Di the 
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planetary radius {Rp = 1.807 ± 0.039 Rjup), and computes 
the Roche shape of the planet which would have the same 
cross-section during the transit as the one we inferred from 
our observations assuming a spherical planet (eclipse model 
of Mandel & Agol 2002, see Section 3.3). 

The results are displayed in Table 5. Rsuh, Rback, Rpoie, 
and Rside are the planetary radii at the sub-stellar point, 
the anti-stellar point, the rotation pole, and on the side, 
respectively. Together, these parameters describe the Roche 
shape of the planet. Raoss is the cross-sectional radius, 
i.e. the radius of the sphere that would have the same 
cross-section as the Roche surface of the planet during the 
transit. It is the planetary radius we derived from our global 
analysis in Section 3.3 (i.e. the observed radius). Rmean is 
the radius of the sphere with the same volume as the Roche 
surface of the planet. The Roche-lobe filling parameter ff 
is defined as Rauh/Rhi, where Rli is the distance of the 
Lagrangian LI point. The asphericity of the planet can 
be quantified by the ratios Rsub/Rside, Rsub/Rpoie, and 
Rside/Rpoie. Finally, the quantity (Rcross/Rmean)® is the 
correction factor that must be applied to the planetary 
density pp derived from our global analysis assuming a 
spherical planet to convert it to the density obtained using 
the Roche model. 

With a Roche-lobe filling parameter of 0.59, 
WASP-121 b is still well away from Roche lobe over¬ 
flow, despite being significantly deformed. It would 
nonetheless be interesting to search for potential signa¬ 
tures of planetary material surrounding WASP-121 (e.g. 
excess transit depths in the near-UV). Such signatures 
have indeed been possibly detected for WASP-12 by 
Fossati et al. (2010), Haswell et al. (2012), and Fossati 
et al. (2013), although Budaj (2011) reported a value 
of only 0.61 for the Roche-lobe filling parameter of 
WASP-12 b. If we use the Rsub/Rpoie ratio to quantify 
the departures from the sphere, we find that WASP-121 b 
(Rsub/Rpoie=1.124) is one of the most distorted planets 
known to date, alongside WASP-12b (7?sub/Rpoie=1.138, 
Budaj 2011), WASP-19b (Rsub/Rpoie=1.124, Budaj 2011), 
and WASP-103b (Raub/Rpoie=1.120, Southworth et al. 
2015). Rmean being more representative of the physical size 
of the planet than Rcross, we adopt it as our revised value 
for the planetary radius (1.865 ± 0.044 Rjup) and include 
it, as well as the subsequent revised value for the planetary 
density (0.183 ± 0.016 pjup), in Table 4. We note that 
Rmean should be used when comparing the radius of the 
planet with theoretical models (e.g. Fortney et al. 2007), 
while Rcroaa should rather be employed when interpreting 
transmission or occultation data. 

Several works showed that hot Jupiters’ radii correlate 
well with their incident irradiation (see e.g. Demory & 
Seager 2011, Enoch et al. 2012, or Weiss et al. 2013). 
Fig. 12 (right) shows the position of WASP-121 b in an 
irradiation-radius diagram for the same sample of transiting 
planets as previously. With a radius of 1.865 ± 0.044 Rjup 
and an irradiation of ~7.1 10® erg s“^cm“®, WASP-121 b 
joins the handful of extremely irradiated planets with 
super-inflated radii. Its radius is actually significantly 
larger than the value of 1.50 ± 0.03 Rjup predicted by the 
equation derived by Weiss et al. (2013) from a sample of 
103 transiting planets with a mass between 150 Mq and 
13 Mjup and that relates planets’ sizes to their masses 


Parameter 

Value 

Radius at the sub-stellar point Rgub [-^Jup] 

2.009 

± 

0.072 

Radius at the anti-stellar point Rback [^Jup] 

1.997 

± 

0.069 

Radius at the rotation pole Rpoie [-^Jup] 

1.787 

± 

0.038 

Radius on the side Rside [-^Jup] 

1.828 

± 

0.041 

Cross-sectional radius Rcross [^jup] 

1.807 

± 

0.039 

Mean radius Rmean [Rjup] 

1.865 

± 

0.044 

Roche-lobe filling parameter ff 

0.591 

± 

0.040 

Rsuh / -^side 

1.099 

± 

0.020 

Rsuh/ Rpole 

1.124 

± 

0.026 

Rside/Rpole 

1.023 

± 

0.003 

Density correction factor (Rcross/Rmean)^ 

0.910 

± 

0.013 


Table 5. Parameters describing the shape of WASP-121 b, ob¬ 
tained using the Roche model of Budaj (2011). The errors are 
the quadratic sum of the errors due to each input parameter (a, 
Mt,IMp, and Rp). 

and irradiations. Several physical mechanisms have been 
proposed to explain the inflated radii of hot Jupiters, such 
as tidal heating (Bodenheimer et al. 2001), deposition of 
kinetic energy into the planets from strong winds driven 
by the large day/night temperature contrasts (Showman 
& Guillot 2002), enhanced atmospheric opacities (Burrows 
et al. 2007), reduced heat transport efficiency by layered 
convection inside the planets (Chabrier & Baraffe 2007), 
or Ohmic heating from currents induced through winds in 
the planetary atmospheres (Batygin & Stevenson 2010). 
As the WASP-121 system is quite young (1.5 ± 1.0 Gyr, 
see Section 3.4), tidal circularization and dissipation might 
have occurred recently enough to contribute to the observed 
inflated radius. 

5.2 Orbital evolution of WASP-121 b 

In Section 3.3, we find a sky-projected spin-orbit angle of 
/3=257.8l5 g deg. This high spin-orbit misalignment favors 
a migration of the planet involving strong dynamical events, 
such as planet-planet scattering (e.g. Rasio & Ford 1996, 
Weidenschilling & Marzari 1996, Moorhead & Adams 2005, 
Chatterjee et al. 2008) and/or Kozai-Lidov oscillations (e.g. 
Kozai 1962, Lidov 1962, Wu & Murray 2003, Fabrycky & 
Tremaine 2007). As shown in Fabrycky & Winn (2009), 
the true spin-orbit angle ij), can be computed from jS, ip 
(the orbital inclination), and i* (the stellar inclination). If 
we assume that the stellar rotation period Rrot,* is either 
0.89 days or 1.13 days and we combine this information with 
our measured values for R* and u* sin i*, we obtain i* of 
9.4 ± 0.6 deg and 11.9 ± 0.7 deg, respectively (see Section 
3.2). These values, together with our measured P and ip, 
yield of 89.6 ± 1.1 deg and 90.1 ± 1.2 deg, respectively, 
where the uncertainties are the quadratic sum of the uncer¬ 
tainties due to each input parameter. Both solutions are in 
excellent agreement and place the planet in a nearly polar 
orbit around its star, which is clearly a remarkable configu¬ 
ration as noted by e.g. Lai 2012 and Rogers & Lin 2013. 

The proximity of WASP-121 b to its Roche limit sug¬ 
gests that its orbital evolution is now dominated by tidal in¬ 
teractions with its host star. To assess the future orbital evo¬ 
lution of the planet, we integrate the tidal evolution equa¬ 
tions of Matsumura et al. (2010) forwards in time, assuming 


MNRAS 000, 1-19 (2014) 





WASP-121 b 15 





Figure 13. Future tidal evolution of the orbital semi-major axis (top), the spin-orbit angle (middle), and the stellar rotation period 
(bottom), assuming the current stellar rotation period P*_o is 0.89 d (left panels) or 1.13 d (right panels). The labeled evolutions 
correspond to Q'^ q= 10®, 10^, 10®, and 10®. The red horizontal dashed lines represent the Roche limit. 


that the planet is on a circular orbit (see Section 3.3) and 
that its rotation is now synchronized with its orbit. Un¬ 
der these reasonable assumptions, the planet’s orbital evo¬ 
lution is expected to depend only on the tidal dissipation 
inside the star (Matsumura et al. 2010). The efficiency of 
this tidal dissipation is parameterized by the stellar tidal 
dissipation factor Q*, with a higher Q* meaning a weaker 
tidal dissipation. We assume here that Q* oc l/|2aj* — 2n|, 
where a;* is the stellar angular velocity and n is the planet’s 
mean motion (n = 2ti/P), as recommended by Matsumura 
et al. (2010) for the cases where the stellar rotation is 
not yet synchronized with the orbit (stellar rotation pe¬ 
riod P* 7 ^ P). Under this assumption, Q* changes as 
Q* = Q'q 12 a;* — 2n|/|2a;*,o — 2no|, where 0 index indicates 
the current values. Fig. 13 shows the obtained evolutions for 
the semi-major axis a (top), the spin-orbit angle ip (middle), 
and the stellar rotation period P* (bottom), assuming differ¬ 
ent values between 10® and 10® for the current stellar tidal 
dissipation factor (5*,o S'lid a current stellar rotation period 
P*,o of 0.89 d (left panels) or 1.13 d (right panels). For all 
cases, the model shows that the planet will continue to ap¬ 
proach its host star until reaching its Roche limit, where it 
will be hnally tidally disrupted. Assuming P*,o=0.89/1.13 d, 
we find remaining lifetimes of 4.20/1.22 Myr, 42.0/12.2 Myr, 
420/122 Myr, and 4.20/1.22 Gyr for Q/ q values of 10®, 10^, 
10®, and 10®, respectively. We note that for Q/,o = 10®, the re¬ 


maining lifetime of the planet would be <1% of the youngest 
possible age of the system (500 Myr, see Section 3.4), thus 
giving a low probability for the planet to be detected now. 
This result could favor (3/,o values > 10^ for WASP-121, 
which would be in line with the expectation that hot stars 
(Teff > 6250 K) would have low tidal dissipation efficien¬ 
cies due to their thin or quasi-nonexistent convective en¬ 
velopes (see e.g. Winn et al. 2010, Albrecht et al. 2012). For 
WASP-121, we have Tea = 6460 ± 140 K and a thin con¬ 
vective envelope starting at 0.88 P*. 

Combining the Matsumura et al.’s expression for 
da/At and the Kepler’s third law, we can calcu¬ 
late the current rate of orbital period change : 
(dP/dt)o = -0.0035 (10®/(3/_o) s yr“^ for P*,o=0.89 d and 
(dP/dt)o = -0.0118 (10®/Q'*’o) syr-i for P*,o=1.13 d. We 
can then estimate how long it would take to observe sig- 
nihcant transit timing variations (TTVs) due to the orbital 
decay of the planet using, e.g., the equation (7) of Birkby 
et al. (2014). We note that this estimation assumes a con¬ 
stant rate of orbital period change dP/dt=(dP/dt)o. Fig. 14 
shows the evolution of the shift in the transit time as a func¬ 
tion of time for the different and for P*^o=0.89 d (top) 

or P*, 0=1.13 d (bottom). Assuming that a timing accuracy 
of ~20 s is achievable with current instrumentation (see 
e.g. Gillon et al. 2009) and that P*,o=0.89/1.13 d, transit 
timing variations could be detected at 3a after ~ll/6 yrs 
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Figure 14. Shift in the transit time of WASP-121 b as a func¬ 
tion of time, computed for different values of Q'^ q, assuming the 
current stellar rotation period is 0.89 d (top) or 1.13 d (bottom). 
The labeled evolutions correspond to ^=10®, 10^, 10®, and 
10^. The red horizontal line represents the 3cr detection limit, 
assuming a timing accuracy of 20 s. 


for = ~35/19 yrs for Q*o=10^, ~110/60 yrs for 

(5To=10®, and ~348/189 yrs for (3i_o=10®, respectively. For 
low values of the decay of WASP-121 b’s orbit could 
thus be detectable on the decade timescale. Alternatively, 
the non-detection of TTVs in the WASP-121 b system after 
a certain amount of time would allow to put a lower limit 
on and help constrain tidal evolution theories. 


5.3 Validity of the spherical approximation for 
the star WASP-121 

The likely stellar rotation periods of 0.89 days and 
1.13 days would imply, together with our measured value 
for i?*, rotation velocities v^, of 82.9 ± 1.7 kms“^ and 

65.3 ± 1.4 kms“^, respectively (see Section 3.2). Such high 

w* might imply a non-negligible deformation of the star 
due to its fast rotation. Using the Roche model from 
Maeder (2009), this deformation can be expressed as 
R*,eq/R*,poi = 1-1- (1/2) (a;*/a;*,K)^, where w* is the ro¬ 
tation velocity; a;*,K the critical break-up velocity; R*,eq 
the equatorial radius; and R.,,poi the polar radius. We es¬ 
timate the critical velocity from the Keplerian velocity 
^GMi,/R^^eq, for which the centrifugal force equals grav¬ 
itational forces at the equator of the star. For WASP-121 
(Table 4), = 0.17 which gives i?*,eq/R*,poi = 1.015. 

The deformation of the star induced by its possible rapid 
rotation would thus not strongly bias our stellar radius mea¬ 
surement (which has an error bar ~2%), even if it was ob¬ 
tained assuming a spherical shape for the star (eclipse model 
of Mandel & Agol 2002, see Section 3.3). 


5.4 First constraints on the atmospheric 
properties of WASP-121 b 


The large radius of WASP-121 b, its extreme irradiation, and 
the brightness of its host star make it an excellent target for 
atmospheric studies via secondary eclipse observations, with 
theoretical expectations for the planet-to-star IR flux ratio 
>0.05% down to ~0.9 fim. 

By combining seven occultation lightcurves obtained 
with TRAPPIST, we detect the emission of the planet in the 
z'-hand at better than ~4cr, the measured occultation depth 
being 603 ± 130 ppm (see Fig. 5). To make sure that the de¬ 
tected signal does not result from a systematic effect present 
in one or several lightcurves, we followed the method pre¬ 
sented in Lendl et al. (2013). We thus divided our set of seven 
occultation lightcurves into subsets containing all possible 
combinations of three to six lightcurves and performed an 
MCMC analysis on each of them, while keeping all the jump 
parameters except the occultation depth fixed to the values 
derived from our global analysis (Table 4). Fig. 15 presents 
histograms of the derived occultation depths. We can see 
how the solutions converge towards our adopted value as we 
use an increasing number of occultation lightcurves. 

Our measured occultation depth can be translated into 
a brightness temperature Tbr of 3553ti7g K. In this cal¬ 
culation, we considered the planet as a blackbody of tem¬ 
perature Tbr and used a Kurucz model spectrum for the 
star (Castelli & Kurucz 2004). We also assumed that all the 
planet’s z'-band flux arises from thermal emission, which is 
reasonable as Lopez-Morales & Seager (2007) showed that 
reflected light contributes 10 to 20 times less than thermal 
emission at these wavelengths for very hot Jupiters, such 
as WASP-121 b. We then defined our measured occultation 
depth as the product of the planet-to-star area ratio and 
the ratio of the TRAPPIST z' bandpass-integrated planet- 
to-star fluxes (e.g. Charbonneau et al. 2005), and adopted 
as brightness temperature the blackbody temperature that 
yielded the best match to our measured occultation depth. 
The uncertainty on Tbr only accounts for the uncertainty on 
the measured occultation depth, which is the main source 
of error. 

We can compare this brightness temperature to the 
equilibrium temperature of the planet, which is given by 
Teq = Teff a/R* /fl [/ (1 — Ab)]^^"^, where / and Ab are the 
reradiation factor and the Bond albedo of the planet, respec¬ 
tively (Lopez-Morales & Seager 2007). The factor / ranges 
from 1/4 to 2/3, where /=l/4 indicates efficient heat redis¬ 
tribution and isotropic reradiation over the whole planet and 
/=2/3 corresponds to instantaneous reradiation of incident 
radiation with no heat redistribution. The brightness tem¬ 
perature is significantly higher (~7 (j) than the maximum 
equilibrium temperature of 30131^ K obtained assuming 
Ab=0 and /=2/3, which could suggest that the planet has 
a low Bond albedo coupled to a poor heat redistribution ef¬ 
ficiency. This would agree with the trend noted by Cowan 
& Agol (2011) that hottest planets are less efficient at redis¬ 
tributing the incident stellar energy than colder planets, the 
explanation proposed for this trend being that the radiative 
timescale is shorter than the advective timescale for hotter 
planets, causing these planets to reradiate the incident stel¬ 
lar energy rather than advecting it through winds. 
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z'-band ocoultation depth [ppm] 

Figure 15. Histograms of the z'-band ocoultation depths derived 
from the MCMC analyses of subsets of three (green), four (blue), 
five (purple), and six (red) lightcurves chosen among the seven 
TRAPPIST ocoultation lightcurves. 

However, we note that the comparison of brightness 
temperature with equilibrium temperature has only limited 
physical meaning here as our z'-band observations probe 
thermal emission from deep layers (P~l bar) of the plan¬ 
etary atmosphere, that can be hotter than the maximum 
equilibrium temperature (e.g. Madhusudhan 2012). Further¬ 
more, the emission spectrum of a hot Jupiter can deviate 
significantly from that of a blackbody. For instance, the flux 
observed in the z'-band can be noticeably increased (resp. 
decreased) by emission (resp. absorption) due to gaseous 
TiO, which is expected to be abundant if the atmosphere has 
a carbon-to-oxygen (C/O) ratio <1 (0-rich atmosphere)''. 

Our z'-band occultation measurement thus provides a 
first observational constraint on the emission spectrum of 
WASP-121 b. Combined with future observations at longer 
wavelengths (from the ground in the near-infrared J-, H-, 
and K- bands, or from space between 1.1 and 1.7 /rm with 
HST/WFC3 or at 3.6 and 4.5 /rm with Spitzer/IRAC), it 
will allow to gain insights on the planet’s dayside chemical 
composition and thermal structure. 
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Table Al. Priors used in the FASTIS analyses: U{a-,b) represents a Uniform prior between a and 6; A/*(/i;cr^) represents a Normal 
distribution with a mean of p and a width of cr^; V{a; Xmin \ Xmax) represents a Power Law distribution with an exponent cx. computed 
between Xmin a.nd Xmax ; 'P 2 {oii]OL 2 '-,XQ\XTnin\^raax) represents a double Power Law distribution with an exponent oli computed 
between Xmin and xq and an exponent a2 computed between xq and Xmax'i and finally S{a]b) represents a Sine distribution between a 
and b. Notes. In FASTIS, planets are considered as non-self-emitting objects, hence no thermal emission. Only the refiected light is 
considered. Here, we used the albedo of the planet as a proxy for the depth of the occultation. The upper limit in the prior corresponds 
to a geometric albedo of 1 and a brightness temperature of 3500K. The TRAPPIST photometry was performed on focus and no 
stellar contamination is detected within the photometric aperture. We defined a prior for the contamination centered on zero with an 
uncertainty of 1% (positive and negative values are allowed) to account for possible variation of the sky background flux between the 
observations or for stellar brightness variation due to non-occulting spots. The EulerCam photometry was performed slightly out of focus, 
which means the target should be blended with the light from a nearby star about 6.8 magnitudes fainter (R-band), resulting in 

a contamination of 0.2%. The prior for the EulerCam contamination was thus centered on this value. 


Parameter 

Planet 

Planet in Binary 

Triple 

Target parameters 

Effective temperature Teff [K] 

V(6460; 140) 

V(6460; 140) 

^(6460; 140) 

Surface gravity \ogg [cm.s“^] 

^(4.2; 0 . 2 ) 

^(4.2; 0 . 2 ) 

V(4.2;0.2) 

Iron abondance [Fe/H] [dex] 

^(0.13; 0.09) 

V(0.13;0.09) 

V(0.13;0.09) 

Distance d [pc] 

P( 2 . 0 ; 10 ; 10000 ) 

P( 2 . 0 ; 10 ; 10000 ) 

P( 2 . 0 ; 10 ; 10000 ) 

Interstellar extinction E(B-V) [mag] 

W( 0 ; 0 . 1 ) 

W( 0 ; 0 . 1 ) 

W( 0 ; 0 . 1 ) 

Planet parameters 

Radius Rp [R,„p] 

W( 0 ; 2 . 2 ) 

W( 0 ; 2 . 2 ) 

- 

Albedo“ Ag 

W{0;2.5) 

W(0;2.5) 

- 

Binary parameters 

Mass of stellar host M 2 [Mq] 

- 

'P2(-1.3;-2.3;0.5;0.1;2) 

P2(-1.3;-2.3;0.5;0.1;2) 

Mass of stellar companion M 3 [Mq] 

- 

- 

P2(-1.3;-2.3;0.5;0.1;2) 

Orbital parameters 

Orbital period P [d] 

V(l.2749281; 5. 10 -'^) 

^(1.2749281; 5.lO-"') 

V(l.2749281; 5. 10 -"^) 

Transit epoch To [BJD - 2450000] 

^(6635.70839; 0.03) 

^(6635.70839; 0.03) 

V(6635.70839; 0.03) 

Orbital inclination i [°] 

5(60; 90) 

5(60; 90) 

5(60; 90) 

TRAPPIST data parameters (xl2) 

Contamination^ 

V( 0 ; 0 . 01 ) 

V( 0 ; 0 . 01 ) 

V( 0 ; 0 . 01 ) 

Out-of-transit flux 

W(0.9; 1.1) 

W(0.9; 1.1) 

W(0.9; 1.1) 

Jitter 

W{ 0 ; 0 . 1 ) 

W( 0 ; 0 . 1 ) 

W( 0 ; 0 . 1 ) 

EulerCam data parameters (x4) 

Contamination^ 

V( 0 . 002 ; 0 . 01 ) 

V( 0 . 002 ; 0 . 01 ) 

V( 0 . 002 ; 0 . 01 ) 

Out-of-transit flux 

W(0.9;l.l) 

W(0.9; 1.1) 

W(0.9; 1.1) 

Jitter 

W{ 0 ; 0 . 1 ) 

W( 0 ; 0 . 1 ) 

W( 0 ; 0 . 1 ) 

Spectral Energy Distribution parameter 

Jitter [mag] 

W( 0 ; 1 ) 

W( 0 ; 1 ) 
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Table A2. Posterior distributions results of the FASTIS analyses. Notes. This corresponds to an occultation depth of 620 di 140 
ppm. This corresponds to an occultation depth of 2470 di 370 ppm. The error does not account for the uncertainty of the stellar 
models. 


Parameter 

Planet 

Planet in Binary 

Triple 

Target parameters 




Effective temperature Tgfr [K] 

6650 ± 60 

6140 ± 150 

6763lt6o 

Surface gravity log^ [cm.s~^] 

4.245 ± 0.006 


Q 'y^/^ + 0.002 

O. / < U_Q 004 

Iron abondance [Fe/H] [dex] 

0.1 ± 0.1 

n Cl +0.01 

U.Ul_o 22 

-o.2+S:“= 

Distance d [pc] 

261 ± 4.6 

334 ± 12 

578iIo 

Interstellar extinction E(B-V) [mag] 

u.uuy_o.oo7 

U.UU i i _o ,006 

0.048tS;S^t 

Planet parameters 




Radius Rp [Rj„p] 

1.76 ± 0.02 

z.±yyo_o,ooo5 

- 

Albedo Ag 

0.55 ± 0.12“ 

1.14 ± O.IS*” 

- 

Binary parameters 




Mass of eclipse host M 2 [Mq] 

- 

0.947 ± 0.007“ 

1 41 + 0 . Olc 

-^•^^-0.03 

Mass of eclipse companion M 3 [Mq] 

- 

- 

n ocy+o.oosc 

U.OU/ _0.008 

Orbital parameters 




Orbital period P [d] 

1.27492477 ± 5.4 10“'^ 

1.2749280 ± 1.7 10“® 

1.27492495 ± 7.4 lO"'^ 

Transit epoch To [BJD - 2450000] 

6635.7077 ± 0.0001 

6635.7077 ± 0.0034 

6635.7077 ± 0.0002 

Orbital inclination i [°] 

89.131S:™ 

89.90 ± 0.09 

89.78iSi® 

Spectral Energy Distribution parameter 




Jitter [mag] 

0.013 ± 0.013 

0.170 ± 0.072 

0.014tS;Sit 
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